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Abstract: We obtain the general cosmological evolution equations for a classically con- 
sistent theory of bimetric gravity. Their analytic solutions are demonstrated to generically 
allow for a cosmic evolution starting out from a matter dominated FLRW universe and 
relaxing towards a de Sitter (anti-de Sitter) phase at late cosmic time. In particular, we 
examine a subclass of models which contain solutions that are able to reproduce the expan- 
sion history of the cosmic concordance model inspite of the nonlinear couplings of the two 
metrics. This is demonstrated explicitly by fitting these models to observational data from 
Type la supernovae, Cosmic Microwave Background and Baryon Acoustic Oscillations. In 
the appendix we comment on the relation to massive gravity. 
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1. Introduction and summary 

Cosmological observations confirm to a high degree of precision that our universe is ho- 
mogeneous and isotropic and can be described by a single metric using Einstein's theory 
of general relativity with the addition of a cosmological constant Aobs- At the same time, 
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one of the corner stones of this theory is the equivalence principle which, taken together 
with quantum field theory, provides one of the most serious unresolved problems of mod- 
ern theoretical physics. Quantum field theory suggests the reality of vacuum energy and 
naturally provides a cosmological constant Avao while the equivalence principle guarantees 
that this energy gravitates with the same strength as any other source. The main problem 
is to reconcile the expected large value of Avac with the observed small value of Aobs , with- 
out recourse to extreme fine tuning. This is one of the main motivations for considering 
theories of modified gravity. 

An intuitive modification of general relativity is massive gravity which generically also 
suffers from a ghost instability [||. Such theories must contain at least two metrics, say, 
Qf^u and /^jy, where ff^u is treated as fixed and non-dynamical. This extra metric is needed 
to construct non-derivative, non-linear mass terms in the action; since \/—g and g^n = 4 



alone are not adequate for the purpose. Massive gravity has received increased attention 
since a class of massive actions, formulated with a flat /^jy, was proposed in |||, ^ and were 
shown to be ghost free at the completely non-linear levelQ. These were soon generalized 
to arbitrary f^i, and proved to remain ghost free |^, The potential relevance of massive 
gravity to the cosmological constant problem can be argued on a heuristic level, based 
on the Yukawa suppression of the massive amplitudes over large distances. This may be 
expected to weaken gravity over such distances and mimic cosmic acceleration or possibly 
screen out a large cosmological constant. "Self-accelerated" solutions of massive gravity 
with a small cosmological constant have been considered, for example, in [0, ^. However, 
the screaning mechanism does not seem to work without fine tuning 

Here, we will consider a generalization of massive gravity, the bimetric gravity, where 
f^i, is promoted to a dynamical variable with its own kinetic term. The idea of studying 
bimetric gravity for cosmological purposes has a long history [10, 11] and such theories 
are known to admit cosmological vacuum solutions [jl2| . For more recent work on bimetric 
theories see |13, |l^, |l^, |T^, 17, 18, 19| and references therein. In these studies the nonlinear 
bimetric interaction potential is often chosen with the only restriction that it reduces to 
the Fierz-Pauli form |2^, ^ at the linearized level. However, these theories also contain 
the Boulware-Deser ghost Q, similar to massive gravity, and hence are not consistent. It 
has been known that internal consistency imposes restrictions on the kinetic structure of 
bi- metric theories [22|. Recently, a class of ghost free bimetric theories were constructed 
in [23]. Here we will focus on the cosmological implications of these theories. Spherically 
symmetric vacuum solutions of these models were considered in [24]. 

Bimetric theories describe a pair of interacting massless and massive spin-2 fields 
12, ^ which are combinations of the two metrics. The presence of the massless mode 



distinguishes bimetric gravity from pure massive gravity. A priori it is not obvious what 
combination of the spin-2 fields should couple to the matter, except that the coupling 
should not reintroduce the ghost. In this paper we consider the bimetric theory of [23] and 
assume that only one of the metrics, Qfj^u, is coupled to matter. Then, although not a mass 
eigenstate, this metric alone will determine the geodesies and the causal structure of the 
spacetime. 

Looking for cosmological solutions we take both of the metrics to be homogeneous 
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and isotropic and obtain the general cosmological equations for the metric g'^jy. Solving 
these in their most general form involves finding the zero locus of a quartic polynominal. 
This guarantees that solutions of the general equations exist and they are straightforward 
to derive. Prom the general equation we conclude some generic properties of the solution 
inherent to the model. First, they generically asymptote towards a de Sitter or anti- 
de Sitter (AdS) spacetime and hence will in general be able to mimic late-time cosmic 
acceleration. Secondly, one branch of the solutions always allows for an expansion history 
which starts out from an ordinary FLRW like universe at high densities, and hence the 
usual early universe considerations still apply for this solution. 

Displaying the full solutions of the quartic in their analytic form is not very illuminating 
for general values of the parameters of the theory. We study instead the analytic details of 
two simpler classes of solutions, where the quartic is reduced to a quadratic. These solutions 
are shown to exhibit the general behaviour expected of the full solution. In particular, we 
identify and study a "minimal" bimetric model in more detail. In this model, the Friedman 
equation is shown to be completely degenerate with general relativity, up to a possible 
rescaling of Gn, the Newton's constant. Hence, observationally it is indistinguishable from 
the usual concordance model on cosmological scales. We demonstrate this explicitly by 
fitting the model to observations, using data from Type la supernovae (SNe la). Cosmic 
Microwave Background (CMB) and Baryon Acoustic Oscillations (BAO). We further study 
neighbouring models of the minimal model and find that cosmological data favours models 
close the minimal one. 

The paper is organized as follows: In section ^ we review ghost free bi-metric gravity 
and present the equations of motion for a homogeneous and isotropic ansatz for both 
metrics. In section ^ we discuss the general properties of cosmological solutions and then 
focuses on two more constrained models that allow for solutions close to those of general 
relativity. We specify the model used for the comparison to observations in section Q and 
summarize the cosmological data considered in section ^. Section ^ contains the results 
of the model fit, which are discussed in section |^. Appendix ^ contains the details for 
the equations of motion in our metric ansatz. In Appendix ^ we discuss the Bianchi 
constraints and comment on the relation to some cosmological solutions recently obtained 
in the massive gravity literature. 

2. The bimetric gravity equations in the cosmological ansatz 

In this section we review the structure of ghost free bimetric gravity and its equations of 
motion. We then write the equations of motion using a cosmological ansatz for the metrics. 

2.1 Review of bimetric gravity action and equations of motion 

The action for the metric Qfj^iy interacting with another spin-2 field through a non- 
derivative potential is determined by the requirement of the absence of the Boulware-Deser 
ghost. The most general action of this type, modulo choice of coupling to matter, has the 
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form H, 



S = - ^ J d^x y^^d^ Rig) - ^ J d^x y^^d^ Rif) (2.1) 

4 

+ m'^M^ 1 d^x V- det 5 ^ /3„ Cni^fg^f) + j d^x^-detg Cmig,^)- 

Here, Pn are free parameters and e„(X) are elementary symmetric polynomials of the 
eigenvalues of the matrix X given explicitly by 

eo(X) = l, ei(X) = [X], e2(X) = i([X]2 - [X^]) , 

e3(X) = i([X]3 - 3[X] [X2] + 2[X3]) , e4(X) = det(X) , (2.2) 

where the square brackets denote the matrix trace. The non-trivial point here is the 
appearance of the square root matrix X = vy^/ which is necessary to avoid the ghost at 
the nonlinear level. 



As demonstrated through a nonlinear ADM analysis in [y, 23 1, this action is ghost 
free and contains 7 propagating degrees of freedom; 2 corresponding to a massless spin-2 
graviton and 5 corresponding to a massive spin-2 field. Both g^i, and f^y are combinations 
of these massless and massive degrees of freedom. 

Since e„ ~ (y^ g~^f)^, the parameterize the order of interactions at the nonlinear 
level. As such they are more convenient to use than other parameterizations. Besides the 
two Planck masses, the action ( |2.lD contains five free parameters Of these, /3o and 
(34 parameterize the cosmological contants of g^jy and /^jy, respectively. One combination 
of the remaining /3's gives the mass of the massive mode, leaving us with two extra free 
parameters to characterize the nonlinear interactions. 

For convenience, and also not to violate the equivalence principle in any drastic way, we 
only couple g^u to matter. This metric determines the geodesies and the causal structure 
of spacetime. Apart from this coupling, the action is invariant under the exchange, 

g^f, /3„^/34-n, Mg^Mf, m^^m^M^/M] (2.3) 

where the last replacement is needed due to our asymmetric parameterization of the mass 
scale in terms of m?. 

Setting /33 = in the action ( ^.l| ) eliminates the highest order interaction term in 
y/g~^f- However, in view of (^]^) we still have a cubic order interaction term in y/ f~^g 
which can in turn be eliminated by setting f3i = 0. In this sense, the choice /3i = /Ss = 
leads to the "minimal" bimetric action, which is the simplest in the class. 

Let us now consider the equations of motion and the Bianchi constraints that follow 



from (2.1). Varying the action with respect to g^,y gives the equations of motion Q 



2 

R^.u - y,uR + ^ ^(-l)"73„ [g,xY^'n)AVg^) + 5.a1("„)^(7^)J = T^T^. , (2.4) 

n=0 



9 
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where y^^^^(\/(7 ^/) are defined below. Similarly, varying with respect to f^u gives, 

,2 3 



R^.u - \f^.R + ^ ^(-l)"/34-n [UxY^i)AVr^) + Uxyin)t.i^/F^)\ = , (2.5) 
* n=0 

where the overbar indicate f^u curvatures and we have introduced the dimensionless ratio 
of Planck masses 

r2 



M] 

w 

9 



(2.6) 



Note that ( |2.5| ) is essentially obtainable from (^^) through the replacements ( |2.3| ) . Finally, 
the matrices Y^^^i^) in (2.4) and ( |2.5[ ) are given by (with square brackets denoting the 
trace) |^] , 



y(o)(x) = i, y(i)(x) = x-i[x], 
y(2)(x) = x2 - x[x] + ii([x]2 - [x2]) , 

y(3)(X) = X^ - X2[X] + iX([X]2 - [X^]) - il([X]3 - 3[X][X2] + 2[X3]) 



(2.7) 



As a consequence of the Bianchi identity and the covariant conservation of T^u^ the 
Qfj^i, equation of motion ( p.4| ) leads to the Bianchi constraint, 

3 



V/^ J^(-l)-/3, g,,Y^^^^{y^) + g.xY>^^^{y^) 



n=Q 



0. 



Similarly, the /^jy equation of motion ( |2.5| ) leads to the Bianchi constraint. 



v^^(-i)"/34_„ kAy("„),(/F^) + /.A^4.(v^) 



n=0 



0. 



(2.8) 



(2.9) 



where the overbar indicates covariant derivatives with respect to the /^i^ metric. Both 
of these Bianchi constraints follow from the invariance of the interaction term under the 
diagonal subgroup of the general coordinate transformations of the two metrics and hence 
are equivalent. From now on we explicitly concentrate on ( |2.8| ). 

Before proceeding further, let us briefly remark on the relation to massive gravity. This 
corresponds to freezing the dynamics of ff^i, and hence loosing the corresponding equation 
of motion ( |2.5| ) . Hence in the bimetric theory f^i, is much more constrained than in massive 
gravity. 

2.2 Equations in the cosmological ansatz 

We look for solutions where both of the metrics exhibit spatial isotropy and homogeneity. 
For simplicity, we also assume both metrics to have the same spatial curvature k = 0, ±1. 
Then, modulo time reparameterizations, the most general form for the metrics is, ^ 

f^^dxt'dx" = -X^{t) dt2 + Y^{t) dx2 



(2.10) 



^Modulo non-perturbative solutions that can exist for certain values of the /3„ parameters feq] (see also 
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where 



df2 



(2.11) 



1 — kr"^ 

Obviously, we have used time reparameterizations to set goo = — 1 so that the Qfj^u metric 
is in the usual FLRW form. No other transformations are available to further reduce the 
number of functions in the metric ansatz. Thus we need to keep three arbitrary functions 
to describe the metrics with the specified symmetries. Below we will see that the Bianchi 
constraint (implied by general covariance) enforces ^ = ^ = thus leaving only two 
free functions to work with. 

Now, we write down the Bianchi constraint and the equations of motion for the above 
parameterization of the metrics. For the ansatz (2.10) the Bianchi constraint (p.^), or 
dH), gives. 



3m 



Y y2 • 

/3i + 2-/32 + ^/33 



Y -hX 



0. 



(2.12) 



One way this can be satisfied is by setting to zero the expression within the square brackets. 
This implies solutions where y oc a, with special values for the constant of proportionality. 
As will be evident from the Friedmann equation ( |2.14 ) below, this leads to the ordinary 
general relativistic equations with a cosmological constant of order m?. Further, the special 
values imply a vanishing mass for the massive spin-2 field (as is discussed further in the Ap- 
pendix One can also check that linear metric perturbations around these backgrounds 
are indistinguishable from general relativity. Hence one concludes that this is effectively a 
decoupled class of solutions that do not modify general relativity on any scale. 

The true dynamical constraint is enforced by the vanishing of the expression within 
the round brackets, 

X^i.^, (2.3) 

a da 

Using this result together with the ansatz ( p.lOj ) , the g^^y equations of motion ( |2.4D lead to 
the modified Friedmann equation (for details see Appendix^), 



, a\ k 
-3( - -3^ 



m 



Y Y"^ Y'^ 

/3o + 3/3i- + 3/32^ + /33^ 



Mf °' 



(2.14) 



and the acceleration equation. 



+ m 



/3o + 2/3 




/y2 yy\ ^ Y^Y 

+ P2 — + 2— +/33^- 
\ aa j a'^a 



(2.15) 



Here we have used T\ = T\ = T\, consistent with the symmetries of the spacetime. These 
obviously reduce to the ordinary Friedmann and acceleration equations of cosmology in the 
limit m? — ?■ 0, although the solutions will not always be well defined in this limit. Similarly, 
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the f^u equations of motion ( |2.5D lead to, 



-''f) - 3^ + 52(/54 + 3^3^ + 3/32^ + /3i^ 



0, 



and 



aa 



YY 



d\ k m 
Y ~Y^^W 



(2.16) 



0. (2.17) 



The first of these, the /-Priedmann equation ( 2.16| ), is in general a cubic equation for Y 
so the system can be solved exactly. Further, as a consequence of already having imposed 
the Bianchi constraint the two equations are not independent. Indeed, the /-acceleration 



equation ( 2.17 ) is obtained by acting with (3 + (y/y)9j)/3 on the /-Friedmann equation. 
Thus we only need to consider ( p. 16 ). 

Let us briefly comment on the source structure. In what follows we will assume a 
perfect fluid source T^y = [p + P)Ufj_Uy + Pg^u so that in the rest frame. 



rpU 



-p , T\ = P (no sum implied) 



(2.18) 



Assuming also an equation of state of the usual form, P{t) = wp{t), the continuity equation, 

3 a / „ 1 a \ , ^ 

tells us that, for w 7^ —1, 

pIp 



Po 



ao 



-3(1+10) 



and H 



3(1 + ^1-) 



(2.20) 



Here po is the present day energy density. Now, in an expanding universe any source with 
w > —1 will get diluted as the scale factor grows. Hence, p < and H > 0. 

3. Viable cosmological solutions 

In this section we will consider cosmological solutions in bimetric gravity. In particular we 
concentrate on parameter values for which the solutions are close to cosmological solutions 
in general relativity and hence are not observationally ruled out. 

3.1 General features of solutions 

In order to make the analysis more transparent, we deflne the dimensionless combinations 



Y 



T = - 

a 



p^ = p/3M^m\ (3.1) 
In terms of these, the two Friedmann equations ( |2.14| ) and ( p. 16 ) can be written as 

(3.2) 



^T3 + /32T2 + /3iT + ^ + P. 



k 
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and 

04 _9 03 01 I 02 k , , 

3Af2 3M2 T ^ M2 m2 m2a2 ' ^ ^ 

Subtracting these two equations to eliminate gives in general a quartic equation for T, 

3 V 3M2 y V y V 3 M2 J 3M2 ^ ^ 

This determines T as a function of p^. The analytic solutions are straightforward to derive 
but lengthy to display. For generic values of parameters there is the further complexity 
that in some cases one has to match different branches of solutions in order to keep T real 
for all p, although a real solution always exists. In order to keep the discussion transparent 
we simply state some generic properties here and then consider a few special cases in more 
detail. 

First, at late times in an expanding universe, p generically approaches a constant 
value Pvac, the vacuum energy density2. Then from ( [3.4| ) it is obvious that at late times T 
approaches a constant value and this holds for generic values of parameters. The Friedmann 
equations (|3.2| ) and ( |3.3D then imply that at late times the universe always asymptotes 
towards a de Sitter or anti-de Sitter geometry. 

On the other hand, in the early time limit of — )■ oo (in practice, >> 0n, 0n/M^ , 
but not necessarily including 0q) the linear term in ( |3.4[ ) dominates and forces a solution 
T — 7- 0. The Friedmann equation ( p.2[ ) then implies an evolution for dominated by p^, 
and a cosmological constant contribution 0q, as in general relativity. Of course, even in 
this limit, the full quartic equation has three more solutions (two being complex), but these 
diverge in the limit p+ — oo, as can be inferred from the explicit form of the solution. Thus, 
the solution that vanishes is physically preferable since we then get back ordinary general 
relativity very early in the expansion history. While a priori this may not be necessary, it 
guarantees that ordinary early universe considerations remain valid. 

As a caveat, note that while the discussion above is fomulated in terms of T and holds 
in general, at the end, it is only the behaviour of that concerns us. Later we will display 
an example of an analytic solution where T diverges at early times while the Friedmann 
equation for is completely equivalent to the usual general relativistic equation. 

To summarize, the equations admit generic solutions that start out from a universe 
described by the ordinary general relativistic Friedmann equation with the cosmological 
constant ~ 0Qm?-, and evolve toward a de Sitter (or AdS) universe, with a cosmological 
constant depending on the parameters of the theory. 

3.2 Models with /Ja = 

In the following we consider primarily theories with 0^ = 0. Although our main intent is to 
simplify the analysis of the solutions we note that this choice is of interest as it corresponds 



to neglecting the highest (cubic) order nonlinear interactions in g ^ f in ( |2.lD . Also, in 
the massive gravity limit, which corresponds to freezing the dynamics of f^^^ only for this 



^Of course, we can always choose pvac ~ and include the vacuum energy contribution entirely in the 
Po parameter, or vice versa. 
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choice the authors in |28] find sphericahy symmetric cosmological solutions that exhibit 
the Vainshtein mechanism in a manner consistent with observations^. 

For the choice /Js = 0, (|3.4|) reduces to a cubic equation for T and one can integrate 
the gfMu— Friedmann equation (|3.2| ) for H. Writing down the exphcit solution is not very 
illuminating, and we simply note that it exists for generic f3i. This solution is fitted to 
data in the following sections. 

Here, however, to illustrate the general discussion of section 3.1, we look for simpler 
analytic solutions. Note that the more restrictive choice /^s = and /34 = 3/?2-/Vf^ converts 
( |3.4| ) into a quadratic equation (the other possibility /Ss = 0, /3i = is discussed below). 
Now, solutions exist only for (3i ^ 0, 

3p. + ^0 - 3P2M-^ ^{p,) 

^ = ^'-'^ 

where for convenience we have defined 



^ = ^(3p, + /3o - 3/32M,-2)2 + 12/32 M,-2 _ (3 Q) 

Choosing the solution with the positive sign ensures that T — )• as p — )• 00, recovering stan- 
dard early cosmology as discussed above. For this choice of sign, the Friedmann equation 
(|3.2| ) becomes, 

H' + ^.=Ci^.+rn^C,-!^ i2C, - 3) ^{p.) + ^g^,^ _ ^p^^^p^) , (3.7) 



where Ci^2 a'l'e constants given by 

Wi (^^' ^ ^^^^^ 



^1 = 1^2 ("3/32 + 2/30/32 - 6-^2 



^2 = U^oPl + /3o^/52 + 1541 - + 9^ ) . (3-8) 



I8/32 V ' " M2 M2 M4 

Although this appears to be a highly nontrivial modification of the general relativistic 
Hubble parameter, in the high energy limit of large p we have, 

while as the energy density dilutes towards a constant value p^aci 

where the constant C3 can be read off from the above expressions. This explicitly demon- 



strates the general arguments of section 3.1 



^Explicitly, find that for Pz 7^ 0, the Vainshtein mechanism screens not only the scalar mode but 
also the tensor modes and hence is ruled out. 
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3.3 Models with (3^ = (3^ = 

Another particularly interesting class of solutions arise when both /3i = /^s = 0, and 



(^2 < /34/3M^. As discussed in section 2A, this corresponds to the minimal bimetric model 
where the interactions are of the lowest order simultaneously for both \fg~^ and \fj~^- 
This keeps only the quadratic order nonlinear interactions in both sectors such that the 



action (2.1) looks like a nonlinear action with a mass potential. 



For these values we find the solution 

^,^3,. + ft-3ftM-^ (3.11) 

AM.-^-3A 

and the Friedmann equation (|3.2D is given by 



a2 /34 - 3/32M2 3M2 ^ 3(/34 - 3/32^2) ' ^ ' 

In general, we can also split the energy density into matter and vacuum contributions, 
P = Pm+Pvac- This equation is completely degenerate with the ordinary general relativistic 
equation with a rescaled Planck mass and a shifted cosmological constant. Thus it can be 
easily fitted to cosmological data. From the observational perspective any discrepancy with 
data must then be looked for by examining the corresponding solutions at smaller scales 
where the homogeneous solutions are not appropriate, e.g. cluster scales. To make this 
point more explicit and also consider neighbouring models with /3i 7^ we proceed to fit 
these models against observational data from supernovae, cosmic microwave background, 
and baryon acoustic oscillations in the upcoming sections. 

Note that this solution corresponds to a case where even though T diverges for large p, 
the equation for H has the general relativistic form. This case evades the general arguments 



of section 3.1. 



4. Parameterization of the solution 



The bimetric action ( p.lj) has a number of degenerate parameters. In this section we discuss 
the choice of parameters to facilitate comparison to data. 

4.1 Parameterization used for model fitting 

In order to make the comparison with standard cosmology more transparent, we define 

^^ = §2^ (4-1) 
where Hq is the present day Hubble parameter. We further define the density parameters, 

1^^—^, = J]^^M2f^T3 + /32T2 + /3iT + ^V (4.2) 

3M2/72 ^2^2 \^ 3 3 y ' V ; 

where as usual 

n = n^{i + z)^ + n^{i + zf + nA + ■■■ (4.3) 
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is given in terms of the fluid components, respectively, for radiation, matter and vacuum 



etc., at redshift z = 0. Dividing through by Hq, the Friedmann equation then assumes the 



form 

e'^ = n + nk{i + zf + nr , (4.4) 

where E'^ = 1 at redshift z = hy definition. 

In the action ( |2.lD , the parameter /3o is a cosmological constant for the metric and 
hence degenerate with the vacuum energy pvac- This can be used to set 

/3o = -3/3i - 3/32 - /33 . (4.5) 

The reason for this choice is simply to have similar conventions to those sometimes used 
in massive gravity (where this choice eliminates the contribution of the mass potential to 
the cosmological constant but only for backgrounds where g = f)- This parameter choice 
factorizes fix as 

Ot = M^^^ll [/33(t2 + T + 1) + 3/32(T + 1) + 3/3i] . (4.6) 

Note that Tq = 1 is a consistent normalization that can be achieved on rescaling /^j^ by 
adjusting Mj. This fixes, Q + $7^ = 1 at the present epoch. The Friedmann equation for 
fnu ( |3.3| ) then allows us to determine, e.g., /34 as, 

ft = -(A+3ft)-3(&-^--jj|). (4.7) 

Finally, note that the parameters /3i , /32 , (^3 are degenerate with the mass scale . In 
particular, we can fix, 

/3i = -1-2/32-/33, (4.8) 

to render the Fierz-Pauli mass of the massive fluctuation independent of the /3's. Then for 
this choice, 

ml, = ^^m^ (4.9) 

is the Fierz-Pauli mass of the massive fluctuation when expanding the metrics around a 
common background for canonically normalized fluctuations. For Mf >> Mg the /^j^ 
sector decouples and becomes the mass for the massive fluctuation of g^jy. 

4.2 Model specific considerations 

We consider a model with zero spatial curvature, k = 0, and /33 = 0. With the above 
choices for /3o ,/3i ,/34, these can now be parameterized as, 

/3o = 3 + 3/32 = 6 - 3a 
/3i = -1-2/32 = -3 + 2a 

/34 = 1-/32 + 3^ = a + 3^ (4.10) 
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where we have defined 



a 



1-/32 



(4.11) 



Note that this is the a parameter of [|9[], discussed further in Appendix |B.1| . We can now 
write the two Friedmann equations (|3.2| ) and (^^) in terms of a as 



= ^^{2- a)M^ + (2a - 3)M^T + (1 - a)M^T'^ 



and 



Ml a 
M2 + 3 



a 



Ml 
M2 



^ 2a 

T H 

3 



1 = 0. 



(4.12) 



(4.13) 



Since our goal is to obtain an analytical expression for E{z), we first substitute the expres- 
sion for E'^ given in ( [4.12| ) into (4.13) to obtain 



« Ml , ,.,,2' 



1^2 



+ 



1 



+ (3- 2a)M,^T 
a - + (« - 2)m2 



^ 2a 

T H 

3 



1 = 0. 



(4.14) 



We then solve this cubic equation for T (being careful to pick the solution corresponding 
to To = 1) and put this back into the Friedmann equation ([4.12|) . The resulting solution 
is fitted to data in the next sections. 

In the case that a = 3/2 (corresponding to /3i = discussed earlier), we obtain a 
particularly simple form for the expansion history, 



where, 



E' 



nf 



n 



1 



nf 



+ nf. 



M'^M? 



(4.15) 



(4.16) 



M2 + M2(2 + M2) ■ 

appears as an effective cosmological constant and also contributes to an effective Planck 
mass. To compare with data, we do not include any vacuum contribution in 0, which, at late 
times, is then entirely given by pm- The actual outcome of this evolution equation depends 
on the relation between the scale Mg and the physical Planck mass Mp (or equivalently, 
the Newton constant G^). In general this will be of the form. 



Mp = QM„, 



(4.17) 



where Q is given by the parameters of the theory. To explicitly determine Q, we need 
localized bimetric solutions that are not well understood at present, so we treat it as a 
parameter ^. Then, in view of (4.2) $7 is related to the physical density parameter by 

f^ = J^phy.Q' (4.18) 



*It is easy to determine Q in the linearized theory which will also exhibit a milder verion of vDVZ 
discontinuity. However, the real Q for non-linear solution will be different due to the Vainshtein effect. 



- 12 - 



and the evolution equation becomes, 



= J^phys - ^f) + . (4.19) 

Comparison to data, gives il^f = 0-''' so that one always has r2phys,oQ^ = 1- The limiting 
values of Q are Q = \ {Mg = Mp) and Q'^{1 — ~ ^ which is equivalent to the 

concordance model. 

As some limits, let us consider, 

M^oo f]f^-i^, (4.20) 

and 

M.^oo ^ nf^-Ml_, (4.21) 

This means that we always have a cosmological constant universe, even as either M or 
M* are extremely large, as long as the other mass ratio is of the appropriate size. As an 
example, let us set Q = 1. Then the first of the above limits, eq. ( 4.2C1| ), corresponds to 



M » 1 ^ m» Hq. (4.22) 



This, as can be seen from eq. ( |4.9| ), implies a very large mass for the spin-2 field measured in 
Hubble units. As is straightforward to verify, given that m/ Hq >> 1 we need Mf ~ l.bMg 
in order for to mimic Vl\ ~ 0.7. It is indeed interesting that this bimetric model can 
mimic the concordance model even for a very large mass for the massive spin-2 mode. It 
is not clear, however, how this will effect smaller scale physics. 
The second limit, eq. ( [4.21| ), corresponds to 

» 1 ^ Mf » Mg . (4.23) 

This will effectively decouple the f^i, field, making it a free field determined by the vacuum 
Einstein equations. This in turn makes the fluctuations of g^i, massive with mass m, which 
again can be seen from (|4.9| ). Now, in order for Q'j^ to mimic Qa ~ 0.7, we need m ~ 2.2Hq. 
We also note that the limit 

M << 1 ^ m« Ho, (4.24) 



can be seen from eq. (4.16) to imply ~ 0, so that the limit of vanishing mass for the 



massive spin-2 field (for non-zero Mf) gives no cosmological contribution. 
5. Data 

In this study, we limit ourselves to purely geometrical tests of the expansion history of the 
universe. That is, tests only involving cosmological distances. We defer possible constraints 
involving smaller scale gravity and structure formation to upcoming work. 
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5.1 Type la supernova data 

As being standardizable candles and thus effective distance indicators, Type la supernovae 
(SNe la) are one of the most direct probes we have of the expansion history of the universe. 
In this paper, we use the Union2 compilation of SNe la. This data set contains SNe la 
from, e.g., the Supernova Legacy Survey, ESSENCE survey and HST observations. After 
selection cuts, the data set amounts to 557 SNe la, spanning a redshift range of < z < 1.4, 
analyzed in a homogeneous fashion using the spectral-template based fit method SALT2. 

5.2 Cosmic Microwave Background and Baryon Acoustic Oscillations 

The position of the first Cosmic Microwave Background (CMB) power-spectrum peak, 
representing the angular scale of the sound horizon at the era of recombination, is given by 

where (iA(-z*) is the comoving angular-diameter distance to recombination while the co- 
moving sound horizon at photon decoupling, rg, is given by 

dz, (5.2) 



H{z) 



which depends upon the speed of sound before recombination, c^. Here we use CMB mea- 
surements from the seven-year Wilkinson Microwave Anisotropy Probe (WMAP) observa- 



tions [30|, in this case the WMAP7.2 results reported at lambda.gsfc.nasa.gov, adopting 
the value Ia = 302.56 zt 0.78. We further assume = 1091.12 exactly (variations within 
the uncertainties about this value do not give significant differences to the results). 

Baryon Acoustic Oscillations (BAO) observations are often compared to theoretical 
models using measurements of the ratio of the sound horizon scale at the drag epoch, 
rs(zrf), to the dilation scale, Dv{z). The drag epoch, ^ 1020, is the epoch at which the 
acoustic oscillations are frozen in. A more model-independent constraint can be achieved 
by multiplying the BAO measurement of rs{zd) / Dy{z) with the CMB measurement I a = 
7rdA{z^)/rs{z^), thus cancelling some of the dependence on the physical size of the sound 



horizon scale [31|. In doing this, we are effectively only left with the assumption that 
the observed inhomogeneities in the large scale distribution of galaxies and in the CMB 
temperature refiects the same (redshifted) physical scale. 

In ||3^ , measurements of the ratio rs{zd)/ Dv{z) at two redshifts, z = 0.2 and z = 
0.35, are reported as r^iza)/ Dv{0.2) = 0.1905 ± 0.0061 and rs(zrf)/Di/(0.35) = 0.1097 ± 
0.0036. Before matching to cosmological models, we need to implement a correction for 
the difference between the sound horizon at the end of the drag epoch, Zd ~ 1020, and the 
sound horizon at last-scattering, z^, ~ 1091, the first being relevant for the BAO and the 
second for the CMB. Here, we use rs{zd)/rs{z^,) = 1.0451 ±0.0158, again using WMAP7.2 
results. A possible caveat is that this ratio was calculated using standard cosmology for the 
evolution between the two redshifts. However, we expect this to be a good approximation 
since the redshift difference is relatively small, and the sound horizon at decoupling and 
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Figure 1: Cosmological constraints for a — 3/2 using supernova distances (SN) and the ratio 
of the observed scales of the baryon acoustic oscillations as imprinted in the cosmic microwave 
background and the large scale galaxy distribution (CMB/BAO). In the left panel, a value of 
= 1 is assumed, in the right panel, — 3. 



drag is mostly governed by the fractional difference between the number of photons and 
baryons. Combining this with ij^ gives the numbers we employ in our cosmology fits, 



Dv{0.2) 

dA{z*) 



17.55 ±0.62, 
10.11 ±0.34. 



(5.3) 



L>v(0.35) 

We take into account the correlation between these measurements using a correlation co- 
efficient of 0.337 calculated in [32|. 



6. Results of observational tests 
6.1 a = 3/2 = 0) 

If we fix the value of = Mf/Mg and allow for a cosmological constant, we can (assuming 
a flat universe) fit for the spin-2 mass parameter m and the matter density We expect 
to get a good fit for M = and dm ~ 0.3 since this corresponds to the concordance 
cosmology with Q\ = 0.7. If we are able to get a good fit also when 1, that is 

with no cosmological constant, depends on the value of M^. From Figure |l], we see that 
in the case of M* = 1 (left panel), this is not possible, while in the case of M* = 3 
(right panel), it is, as expected from Eq. ( |4.15| ). In this and all figures hereafter, shaded 
contours shows constraints for SN and CMB/BAO data, respectively, corresponding to 
95 % confidence interval for two parameters. Combined constraints are shown with solid 
lines corresponding to 95 % and 99.9 % confidence intervals for two parameters. 

In what follows, we will set = 1 and i}\ = in order to investigate whether bi- 
metric gravity models can explain the apparent accelerated expansion seen in cosmological 
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Figure 2: Cosmological constraints in the [M, M*] -plane assuming a = 3/2 and VL^ = 1- 



geometrical data. Fitting both and M, we expect to be able to obtain good fits to the 
data as long as the combination VLf = M'^M^/{M'^ + M^{2 + M^)) is close enough to the 
concordance value for the cosmological constant. As discussed previously, this can always 
be achieved if either one of the values is large enough, see Figure ||. 

6.2 General a / 0) 

We now proceed to fit also the value of a (assuming Q\ = 0). For = 3, we obtain the 
result depicted in the left panel of Figure ^ showing a good fit to the data for a ~ 1.4 and 
M ~ 3.0. 

Even more generally, we want to fit a, M and simultaneously. This requires 
care when projecting results on two dimensional surfaces. Since we are able to obtain 
good fits to the data for both M — cc and — oo, we will not be able to cover 
the entire non-negligible probability function in our grid of tested parameter values which 
in principle is required to perform a proper marginalization. This will mostly affect the 
projected constraints for M and where generally a larger parameter space in one of the 
parameters will give more weight to the likelihood function of the other parameter at lower 
values (see Figure ^). In the following, we present results for M™^^ = M™^^ = 10, the 
equivalent of putting a flat prior on the values of M and M* to be in the interval [0, 10]. 
Results in the [a, M] -plane after marginalizing over M* are shown in the right panel of 
Figure Comparing to the left panel, we can see that allowing to vary slightly shifts 
and widens the allowed values of a and M as compared to the case of a fixed value of M^,. 
Marginalizing down to one parameter surfaces (assuming flat prior probablilities) , our data 
constrains the parameter values of the bimetric gravity model to be (at 95 % confidence 
level for one parameter) 

l.l<a<1.5, 2<M<3.5, 1.5<Af, <3. (6.1) 



-16- 




1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 2 4 5 8 10 

M = m/Ho M = m/Ho 



Figure 3: Left panel: Cosmological constraints in the [a, M]-plane assuming 51„i = 1 and Af* = 3. 
Right panel: Constraints in the [a, Af]-plane for rim ~ 1 after marginahzing over Af* with a fiat 
prior in the interval [0, 10]. 

6.3 Interpretation of results 

As expected from the form of the cosmological solutions, the bimetric gravity model is 
perfectly capable of matching geometric cosmological data. For the case of a = 3/2 (/3i = 
0), this was obviously going to be possible for some values, but we find that data clearly 
favours the regions M ~ 2.5,M* > 1.5 or M > 2.5, ~ 1.5 (assuming Mp = Mg 
for simplicity, the more general case having been discussed in section |4.2| ). Allowing for 
an arbitrary a is equivalent to switching on the /3i interactions, with the only restriction 
that we focus here on bigravity models for which /3i = —1 — 2/32, which include massive 
gravity limits of bigravity. In this case, we find that data favours a rather narrow region 
of parameter values. In particular, large deviations from a = 3/2 as well as m = 2>Hq and 
Mf = 2Mg are being disfavoured by data. 

7. Conclusions 

We have investigated cosmological solutions in the unique classically consistent theory of 
bimetric gravity. Under the assumption that both metrics respect the symmetries of spatial 
isotropy and homogeneity, we derived the most general cosmological evolution equations for 
this theory. The generic solution was demonstrated to always allow for a cosmic evolution 
starting out from an ordinary FLRW matter dominated universe while evolving towards a 
de Sitter (or AdS) geometry at late times in the expansion history. We explored further 
a particular class of solutions, corresponding to neglecting cubic nonlinear interactions in 
the potential for the two metrics. A subclass of these solutions, the minimal model, was 
shown to be completely degenerate with the evolution of a universe described by the usual 
cosmic concordance model in general relativity. Using recent data from SNe la, CMB 
and BAO we demonstrated that data favoured regions in parameter space close to this 
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minimal model. Interestingly, assuming Mp ~ Mg, we found that even for a very large 
mass for the massive spin-2 mode can this theory match observations as good as general 
relativity, although under the assumption of order unity values for the relative couplings, 
data favoured a mass of the order of the Hubble scale ra ~ Hq. 

Since the bimetric gravity theory is a highly nontrivial but consistent modification of 
general relativity, it is important to explore its consequences further. In particular, explor- 
ing whether more general classes of solutions than the subclass studied in this paper are 
able to match observations. Since we have demonstrated the existence of bimetric solutions 
that match general relativity on cosmic scales it is important to study perturbations on 
these solutions and also obtain the corresponding small scale solutions in order to see what 
this theory predicts for observations on e.g. cluster and galaxy scale. Another important 
issue is the determination of the observed Planck scale in terms of the parameters of the 
theory. 
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Note Added: When the writing of this paper was being finalized, a paper with similar 
intent appeared |27] with some overlap with section 3 of the present paper. After this paper 
was submitted to the arXiv, appeared which also considers cosmological implications 
of bimetric gravity 



A. Derivation of the equations of motion 



In this appendix we give some details on the derivation of the equations of motion, equations 
( 2.14 ), ( 2.151) , ( 2.16 ), and (|2.17] ). We insert our ansatz for the two metrics ( 2.10|) , into the g- 
and /-equations of motion ( ^.4| ) and (p.5|), respectively. With our ansatz the 00-component 
of the g^y equation of motion ( ^.4[) becomes 



.a 



k 



and the 11-component reads 



2 K-^ + m' [l3o + 



Y 



Y' 



y3 



-3^ - 3^ + /3o + 3/3i- + 3/32^ + /33^ = 



7^0 _ 



Y 

2— + X 
a 



T2 XY 
^ + 2 



XY^ 

P2 + -^P3 



Now we insert the solution of the Bianchi constraint, 

X = ^ = —, 
a da ' 



(A.l) 



Mf 



(A.2) 



(A.3) 
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to rewrite (|A.2| ) as 



-2 2 - ^ + /^o + 

In order to compute the contribution from the kinetic terms to the f^y equations of mo- 
tion ( p.5[ ), we need the expression for the Ricci tensor R^iy{f) in terms of the two scale 
factors. The non-vanishing Christoffel symbols for /^j^ are 

_ X _ YY _ YYr^ _ YYr^sm^O 

ioo-^, ^ii-X2(l-A;r2)' ^22-^^, i 33 - ^2 ' 

r - — - - -rn - kr'^) r'^ - - i 

j-oi-y' ii~i_ytr2' 22- ^^'^J) -Li2-J^i3-^' 

= -rsin2 6'(l - fcr^) , rfg = - sin 6* cos , T^g = cot 6* . (A.5) 



Y Y 


/3i + 


y2 




2- + - 




+ 2— 


a a 


a2 


aa 



YY2 > 
/32 + ^-rh 



Ml 



(A.4) 



From these we compute the non-vanishing components of the Ricci tensor, 

1 / YY X \ 



^22 = - '—jT + + ^^^'y 



Thus, the curvature scalar is given by 



6 / y Yx y2 /.j^2 \ 

^=^(^7-7^ + 72+l^j- (A.7) 
The Einstein tensor then has the nonvanishing components 
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For the equations of motion, we raise one index with f^'^ and consider the 00- as well as 
the 11-component . The 00-component of the equations of motion is then obtained as 



, Y'^ km? 



3 2 

^/3i+ 3^/32 + 3^/33 + /34 



(A.9) 



whereas the 11-component reads 











2 




m 




Ml 




into O 



k 



Yx y2\ 

'y^ + y2 - y2 



a \ ^ 2a 

Pi + I y2 + 



/32+(2^ + l)/33 + /34 



(A.IO) 



_ qA _l ZA 
y2 -^^2 + M2 



3 2 

+ 3^/32 + 3^/33 + /34 



ad a? 



k m? 

''^yy~y^~y^^mI 



2aa\ 



a a 



(A.ll) 



(A.12) 



We now observe that acting with (3 + (Y/Y)dt)/3 on ( |AlI| ) gives ( [Al2D . Thus, the two 
equations are equivalent. 



B. Comments related to massive gravity 

In order to fascilitate a comparison with other works on both bimetric and massive gravity, 
using the consistent interaction term of the action ( ^.iD , we note here the relation between 
the different parameters of the action that has been established in the recent literature on 
massive gravity, in particular |P, We also study the nature of the Bianchi constraint 
and clarify some points related to cosmological solutions recently obtained in the massive 
gravity literature. 

B.l Parameters of massive gravity 

In massive gravity, the metric that does not couple to matter is regarded as a fixed back- 
ground, usually taken to be a flat reference metric (although this is not necessary). This 
implies removing the kinetic strength Mf. One must also fix /3o = — 3/3i — 3/32 — 1^3 in 
order to cancel terms linear in perturbations around this background. In order for m? to 
correspond to the mass of the massive spin-2 mode one has to fix also /3i = — 1 — 2/32 — f^s- 
This effectively eliminates three parameters. Since /34 does not enter the equations of mo- 
tion for (jf^jy one can consistently eliminate this parameter as well. Thus, out of the five /3„ 
parameters of the general bimetric theory only two are important for massive gravity. One 
is left with four free parameters (including the Planck mass and the spin-2 mass). The 
remaining two parameters in the interaction can be chosen arbitrarily and (at least) two 
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conventions have occured in the hterature. First, we note the relations between /3„ and the 
On used in S 



/3o = 6 - 4a3 + a4 , 
/3i = —3 + 803 — 0:4 
= I - 20.3 + (54 , 

/^s = as - 04 . 



(B.l) 



The relation between these and the parameters of [0] is 



"3 = 
0:4 = 

B.2 The Bianchi constraint 

We recall the Bianchi constraint ( |2.12| ) 



-3a3 = 6C3 , 
12a4 = -48d5 



(B.2) 



3m 



Y 



Y2 



/3i + 2-/32 + ^/33 



0. 



(B.3) 



Clearly we can enforce this by looking for solutions for Y which force the left bracket to 
vanish. This will however only result in the ordinary general relativistic equations for g^^, 
with the addition of a cosmological constant proportional to as is evident from ( 2.14| ), 
independent of any dynamics for ff^i, . As such they represent a particular class of screening 
solutions with special values for the parameters of the action. 

The reason these values are special can be quantified further. If we consider perturba- 
tions around solutions where the metrics are proportional, i.e. 



we have that 



g 



-1 



g = 9 + f = Cg + 6f , 



f^C + -g-\5f-C5g)^C + 5M, 



(B.4) 



(B.5) 



where gSM defines the massive fiuctuations up to a constant of proportionality. Using this 



and expanding the interaction term in the equations of motion (2.4) to linear order we 
obtain (excluding linear contributions from the cosmological term) 



(A + 2CP2 + C2/33) [g (Tr(<5M) - 5M)], . . (B.6) 



n=0 



This is just the Fierz-Pauli mass contribution to the linearized equations of motion. Com- 
paring this to the Bianchi constraint ( |B.3| ) we see that forcing the left bracket of the Bianchi 
constraint to vanish for proportional metrics is equivalent to choosing a constant of pro- 
portionality such that the Fierz-Pauli mass term for the fluctuations vanish. A similar 
conclusion holds also for the cosmological solutions when the spatial metrics are propor- 
tional but this requires more work to show than the simple example discussed here. 
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For our ansatz ( 2.10 ), the Bianchi constraint (B^) allows to be less strict and only 
demand that the spatial part of the metrics are proportional, i.e. Y{t) = Ca{t). Then by 
forcing the left bracket of the Bianchi constraint to vanish we can find bimetric solutions 
that read 



,2 2... sa^AiV 



where X is determined from the ff_iu equations of motion (^j^) and the constant C is 
determined from the Bianchi constraint ([B.3|) to be given by 



Although the existence of these solutions might appear interesting, in the bimetric theory 
where we include dynamics for f^u, they do not modify general relativity in any respect 
apart from the contribution of a cosmological source given explicitly by the addition 



/5o + ^ - 3/3i/32 ± 2^\//5i - A/Ss T 2/3 J (31 - M3) , (B.9) 

P3 V P3 PS ^ J 



to the usual Friedmann and acceleration equations of general relativity without a cosmo- 
logical constant. From our previous reasoning this can be understood as a consequence of 
using the parameters of the general theory to impose the vanishing of the massive fluctua- 
tions. Hence, the theory contains only massless spin-2 fluctuations and must be equivalent 
to ordinary general relativity. Since we keep all parameters arbitrary such solutions could 
be expected, and a very similar conclusion was reached when examining spherically sym- 
metric solutions in the consistent bimetric theory of gravity in [p^]. 

If we do not include a kinetic term for as in the massive gravity setups, these 
solutions do allow for regular cosmological evolution for g^^, but again only contribute 
with a constant source addition (c.f. (^)). For such a scenario X is not given as in (|B.7| ), 
since we have used the f^i, equations of motion to derive that solution. Indeed, without 
the dynamical term for f^^, X can be an arbitrary function of time and in particular it 
is possible to choose X such that ff^i, represents an open chart of the Minkowski metric, 
as was recently demonstrated in Q. From this perspective it is also clear why it is not 
possible to find solutions with positive or zero spatial curvature for a non-dynamical f^u, 
there is simply no representation of Minkowski space as a homogeneous and isotropic space 
with nontrivial scale factor for these cases. 

B.3 X = constant 

More generally, we recall the true dynamical Bianchi constraint ( ^.13 ), 



a da 



obtained from the vanishing of the right bracket of (|B.3| ). Note that this encodes also the 
Y oc a solutions obtainable from enforcing the vanishing of the left bracket, but in this case 
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will enforce /^j^ oc 5^,^. Thus, for arbitrary parameters of the theory ( B.IO ) contains all 
non-trivial information about the nature of the constraint. This allow for the simple class 
of solutions 

Y = Cia + C2 T = Ci + — . (B.ll) 

a 

where Ci^2 are constants. The Friedmann equation for g^i, ( |2.4D is then given by 

+ (/32 + /33Ci)4 + |^. (B.12) 

Qj (j Qi 

For C2 = this clearly only amounts to a cosmological constant contribution to the usual 
general relativistic equations. Moreover, if Ci = 1 we see that the choice /3q = — 3/3i — 



3/32 — /^s completely eliminates this contribution (c.f. discussion in Appendix B.l ). 

More interestingly, if C2 7^ the addition to the general relativistic Friedmann equation 
is precisely of the form to add extra fluid components of the conventional type. From the 
perspective of the usual formulation of massive gravity then, where there is no dynamical 
equation for f^^, these solutions can be interesting since they can exist for an arbitrary 
source and give a contribution to all the fluid components except for radiation. Solutions 
of this type were recently discussed in that context in |34]. The added degeneracy in 
parameter space is however a rather unpleasant feature. 

In the context of bimetric theory that we are considering here, the equations of motion 
for /^jy ( |3.3D imply that these solutions will constrain the source, as is evident from the 



quartic equation (3^). In fact, they will tell us that the source will contain terms of the 
form 

(B.13) 



Cia + C2' Cia + C2' a{Cia + C2) ' a^{Cia + C2)' 
up to various multiplicative constants. For C2 = 0, they are of the standard matter fluid 
type but then only the constant contribution can remain, as we have also remarked on 



earlier. For Ci = (such that f^i, is flat) the Bianchi constraint (B.IO) tell us that either 
X = or a is a constant and that there can be no cosmological evolution and hence this 
case is not interesting for our purposes. This conclusion was reached also in the massive 
gravity context in |^5[. In that case however no dynamics where considered for f^^, and 
the possibility of finding solutions by enforcing the non-dynamical Bianchi constraint was 
missed due to a too restrictive parameterization of Minkowski space. 

If we do proceed and eliminate the source in the above prescribed manner we end up 
with an evolution equation where the Hubble expansion is driven by terms of exactly the 
same form as in ( |B.13 ). Although it is intriguing that these terms actually do scale in 



accordance with the usual matter fluid components for large values of o, for the purposes 
of this paper we do not want to impose any such restrictions on the matter sector. 

References 

[1] D. G. Boulware and S. Deser, Phys. Rev. D 6 (1972) 3368. 



-23- 



C. de Rham and G. Gabadadze, Phys. Rev. D 82 (2010) 044020 [arXiv: 1007.0443 [hep-th]]. 

C. de Rham, G. Gabadadze and A. J. Tolley, Phys. Rev. Lett. 106 (2011) 231101 

[arXiv:1011.1232 [hcp-th]]. 

S. F. Hassan and R. A. Rosen, arXiv:1106.3344 [hep-th]. 

S. F. Hassan, R. A. Rosen and A. Schmidt-May, arXiv:1109.3230 [hep-th]. 

S. F. Hassan and R. A. Rosen, arXiv:1111.2070 [hep-th]. 

C. dc Rham, G. Gabadadze, L. Heisenberg and D. Pirtskhalava, Phys. Rev. D 83 (2011) 
103516 [arXiv:1010.1780 [hep-th]]. 

A. E. Gumrukcuoglu, C. Lin and S. Mukohyama, JCAP 1111 (2011) 030 [arXiv: 1109.3845 

[hep-th]]. 

S. F. Hassan and R. A. Rosen, JHEP 1107 (2011) 009 [arXiv:1103.6055 [hep-th]]. 
N. Rosen, Phys. Rev. 57 (1940) 150. 

N. Rosen, In *Erice 1975, Proceedings, Topics In Theoretical and Experimental Gravitation 

Physics*, New York 1977, 273-294 

A. Salam and J. A. Strathdee, Phys. Rev. D 16 (1977) 2668. 

N. Arkani-Hamed, H. Georgi and M. D. Schwartz, Annals Phys. 305 (2003) 96 

[hcp-th/0210184]. 

M. Banados, A. Gomberoff, D. C. Rodrigues and C. Skordis, Phys. Rev. D 79 (2009) 063515 
[arXiv:081 1.1270 [gr-qc]]. 

Z. Berezhiani, D. ComelU, F. Nesti and L. Pilo, Phys. Rev. Lett. 99 (2007) 131101 
[hep-th/0703264 [HEP-TH]]. 

D. Bias, D. Comelli, F. Nesti and L. Pilo, Phys. Rev. D 80 (2009) 044025 [arXiv:0905.1699 
[hep-th]]. 

N. Boulangcr, T. Damour, L. Gualtieri and M. Henneaux, Nucl. Phys. B 597 (2001) 127 
[hep-th/0007220]. 

T. Damour and I. I. Kogan, Phys. Rev. D 66 (2002) 104024 [hep-th/0206042]. 
M. Milgrom, Phys. Rev. D 80 (2009) 123536 [arXiv:0912.0790 [gr-qc]]. 
M. Fierz and W. Pauli, Proc. Roy. Soc. Lond. A 173 (1939) 211. 
W. PauU and M. Fierz, Helv. Phys. Acta 12 (1939) 297. 

C. Aragone and S. Deser, Nuovo Cim. A 3 (1971) 709. 
S. F. Hassan and R. A. Rosen, arXiv:1109.3515 [hep-th]. 

D. Comelli, M. Crisostomi, F. Nesti and L. Pilo, arXiv:1110.4967 [hcp-th]. 

C. J. Isham, A. Salam and J. A. Strathdee, Phys. Rev. D 3 (1971) 867. 

T. Damour, I. I. Kogan and A. Papazoglou, Phys. Rev. D 66 (2002) 104025 
[hep-th/0206044]. 



[27] M. S. Volkov, arXiv: 11 10.6153 [hep-th]. 



- 24 - 



[28] K. Koyama, G. Niz and G. Tasinato, Phys. Rev. D 84 (2011) 064033 [arXiv: 1104.2143 
[hep-th]]. 

[29] R. Amanullah, C. Lidman, D. Rubin, G. Aldering, P. Astier, K. Barbary, M. S. Burns and 
A. Conley et al, Astrophys. J. 716 (2010) 712 [arXiv:1004.1711 [astro-ph.CO]]. 

[30] E. Komatsu et al. [WMAP Collaboration], Astrophys. J. Suppl. 192 (2011) 18 
[arXiv:1001.4538 [astro-ph.CO]]. 

[31] J. SoUerman, E. Mortsell, T. M. Davis, M. Blomqvist, B. Bassett, A. C. Becker, D. Cinabro 
and A. V. Filippenko et al., Astrophys. J. 703 (2009) 1374 [arXiv:0908.4276 [astro-ph.CO]]. 

[32] B. A. Rcid et al. [SDSS Collaboration], Mon. Not. Roy. Astron. Soc. 401 (2010) 2148 
[arXiv:0907.1660 [astro-ph.CO]]. 

[33] D. Comelli, M. Crisostomi, F. Nesti and L. Pilo, arXiv:1111.1983 [hep-th]. 

[34] A. H. Chamseddine and M. S. Volkov, Phys. Lett. B 704 (2011) 652 [arXiv: 1107.5504 

[hep-th]]. 

[35] G. D'Amico, C. de Rham, S. Dubovsky, G. Gabadadze, D. Pirtskhalava and A. J. ToUey, 
arXiv:1108.5231 [hep-th]. 



-25- 



